#include <stdio.h>
#include <stdlib.h>

/*
** The following program is for creating the 21 features from the sequence data
*/
const int MAX_LINE_LEN=610;
//const int LINES_IN_INPUT_FILE=313;
const int LINES_IN_INPUT_FILE=385;

//char infilename[] = "Dubchak_NN_fasta_train_seq.txt";
char infilename[] = "Dubchak_NN_fasta_test_seq.txt";
char temp_filename[] = "temp.txt";
char outfilename[] = "hydro_new_385fcm.txt";
//char outfilename[] = "hydro_new6_test.txt";
//char outfilename[] = "hydro_385_old.txt";
char outfilecomposition[] = "composition_new_385fcm.txt";

void transform_file();
void pure_composition();

int main()
{
	int i,j,count,k,m,linenum;
	char ch;
	
	int counter_composition[3];
	
	int transition12;
	int transition23;
	int transition31;

	float *distribution;
	distribution = (float*) malloc(sizeof(float)*15);
	if(distribution==NULL)
	{
		printf("Memory allocation problem at distribution\n");
		exit(1);
	}
	
	int *line;	
	line = (int*) malloc(sizeof(int)*MAX_LINE_LEN);
	if(line==NULL)
	{
		printf("Memory allocation problem at line\n");
		exit(1);
	}
	
	pure_composition();
	transform_file();

	FILE *fp_in = fopen(temp_filename,"r");
	FILE *fp_out = fopen(outfilename,"w");
	if((fp_in==NULL)||(fp_out==NULL))
	{
		printf("Error opening file");
		exit(1);
	}
	
	// for each line do the following computation and write into a file
	for(linenum=0; linenum<LINES_IN_INPUT_FILE; linenum++)
	{
		// initialize all arrays first
		for(i=0; i<MAX_LINE_LEN; i++)
			line[i]=0;
		for(i=0; i<15; i++)
			distribution[i]=0.0;
		for(i=0; i<3; i++)
			counter_composition[i]=0;
		transition12=transition31=transition23=0;
	
		//fill line array
		j=0;
		while( (ch=fgetc(fp_in)) !='\n')
		{
			line[j] = ch-'0';
			j++;
		}
		// line had j chars
		
		//fill counter_composition array
		for(i=0; i<j; i++)
		{
			if(line[i]==1)
				counter_composition[0]++;
			else if(line[i]==2)
				counter_composition[1]++;
			else if(line[i]==3)
				counter_composition[2]++;
		}

		//count is the length of current line (actually it is equal to j right now)
		count = counter_composition[0]+counter_composition[1]+counter_composition[2];

		//compute transition values
		for(i=0; i<(count-1); i++)
		{
			if(((line[i]==1)&&(line[i+1]==2))||((line[i]==2)&&(line[i+1]==1)))
				transition12++;
			else if(((line[i]==1)&&(line[i+1]==3))||((line[i]==3)&&(line[i+1]==1)))
				transition31++;
			else if(((line[i]==2)&&(line[i+1]==3))||((line[i]==3)&&(line[i+1]==2)))
				transition23++;
		}

		
		// now compute 15 features
		// for each of the groups
		for(i=0;i<3;i++)
		{
			j=0;  //elements of this group seen so far
			
			//first occurence of this group...
			while(line[j]!=(i+1))
				j++;
			distribution[5*i] = (float)(j+1)*100/count;

			//first 25% occurence of this group...
			k = counter_composition[i]/4;
			//read the line till we read first k elems of this group.
			for(m=0,j=0; (m<count)&&(j<k); m++)
			{
				if(line[m]==(i+1))
					j++;
			}
			distribution[5*i+1] = (float)m*100/count;

			//first 50% occurence of this group...
			k = counter_composition[i]/2;
			//read the line till we read first k elems of this group.
			for(m=0,j=0; (m<count)&&(j<k); m++)
			{
				if(line[m]==(i+1))
					j++;
			}
			distribution[5*i+2] = (float)m*100/count;

			//first 75% occurence of this group...
			k = 3*counter_composition[i]/4;
			//read the line till we read first k elems of this group.
			for(m=0,j=0; (m<count)&&(j<k); m++)
			{
				if(line[m]==(i+1))
					j++;
			}
			distribution[5*i+3] = (float)m*100/count;

			// 100% occurence of this group...
			k=counter_composition[i];
			//read the line till we read first k elems of this group.
			for(m=0,j=0; (m<count)&&(j<k); m++)
			{
				if(line[m]==(i+1))
					j++;
			}
			distribution[5*i+4] = (float)m*100/count;

		}

		// now write all this information on one line of output file
		//first composition, then transition, finally ditribution
		float temp;
		for(j=0; j<3; j++)
		{
			temp = (float)counter_composition[j]*100/count;
			//printf("  %4.1f",temp);
			fprintf(fp_out,"  %4.1f",temp);
		}
		
		temp = (float)transition12*100/(count-1);
		fprintf(fp_out,"  %4.1f",temp);
		temp = (float)transition31*100/(count-1);
		fprintf(fp_out,"  %4.1f",temp);
		temp = (float)transition23*100/(count-1);
		fprintf(fp_out,"  %4.1f",temp);
		
		for(i=0; i<15; i++)
		{
			//if((i==5)||(i==10))
			//	fprintf(fp_out," ");
			fprintf(fp_out," %5.1f",distribution[i]);
		}

		//the follwing line is required only for comparing results with dubchuk et. al.
		//fprintf(fp_out,"  %d",count);
		
		fprintf(fp_out,"\n");
		
	} // do all this for next input line

	fclose(fp_in);
	fclose(fp_out);

	free(line);
	free(distribution);

/*	//fp_in is the temp file. remove this file.
	fp_in = fopen(temp_filename,"w");
	fclose(fp_in);
	*/
} 


/*
** The following function is for converting the sequence file to groups ****
*/
void transform_file()
{
	
	FILE *fp_in = fopen(infilename,"r");
	FILE *fp_out = fopen(temp_filename,"w");
	char ch;

	int first=0;

	while((ch=fgetc(fp_in))!=EOF)
	{
		//the line containing > sign is copied as it is
		if(ch=='>')
		{
		//	fputc('>',fp_out);
			while((ch=fgetc(fp_in))!='\n')
				;
		//		fputc(ch,fp_out);
			if(first==0)
				first=1;
			else
				fputc('\n',fp_out);
		}
	
	/*
		//for the actual sequence, we put the group nos (due to avg of 23 scales)
		if((ch=='K')||(ch=='E')||(ch=='D')||(ch=='Q')||(ch=='N'))
			fputc('1',fp_out);
		else if((ch=='R')||(ch=='T')||(ch=='W')||(ch=='G')||(ch=='H')||(ch=='A')||(ch=='S'))
			fputc('2',fp_out);
		else if((ch=='Y')||(ch=='V')||(ch=='C')||(ch=='M')||(ch=='P')||(ch=='L')||(ch=='F')||(ch=='I'))
			fputc('3',fp_out);
	*/
/*
		// for the avg of 6 scales
		if((ch=='S')||(ch=='Q')||(ch=='D')||(ch=='R')||(ch=='E')||(ch=='K'))
			fputc('1',fp_out);
		else if((ch=='A')||(ch=='V')||(ch=='H')||(ch=='G')||(ch=='W')||(ch=='T')||(ch=='N'))
			fputc('2',fp_out);
		else if((ch=='F')||(ch=='I')||(ch=='L')||(ch=='P')||(ch=='M')||(ch=='C')||(ch=='Y'))
			fputc('3',fp_out);
*/
/*
		// for the kyte dolittle scale
		if((ch=='N')||(ch=='Q')||(ch=='D')||(ch=='R')||(ch=='E')||(ch=='H')||(ch=='K'))
			fputc('1',fp_out);
		else if((ch=='P')||(ch=='S')||(ch=='T')||(ch=='G')||(ch=='W')||(ch=='Y'))
			fputc('2',fp_out);
		else if((ch=='A')||(ch=='I')||(ch=='L')||(ch=='F')||(ch=='M')||(ch=='C')||(ch=='V'))
			fputc('3',fp_out);
*/  

		// for the fcm scale
		if((ch=='D')||(ch=='E')||(ch=='K')||(ch=='Q')||(ch=='R'))
			fputc('1',fp_out);
		else if((ch=='A')||(ch=='G')||(ch=='H')||(ch=='N')||(ch=='S')||(ch=='T')||(ch=='V')||(ch=='W'))
			fputc('2',fp_out);
		else if((ch=='P')||(ch=='I')||(ch=='L')||(ch=='F')||(ch=='M')||(ch=='C')||(ch=='Y'))
			fputc('3',fp_out);
/*

	// for testing the code: the follwoing lines will generate the same groups as dubchuk...
		if((ch=='R')||(ch=='K')||(ch=='E')||(ch=='D')||(ch=='Q')||(ch=='N'))
			fputc('1',fp_out);
		else if((ch=='G')||(ch=='A')||(ch=='S')||(ch=='T')||(ch=='P')||(ch=='H')||(ch=='Y'))
			fputc('2',fp_out);
		else if((ch=='C')||(ch=='V')||(ch=='L')||(ch=='I')||(ch=='M')||(ch=='F')||(ch=='W'))
			fputc('3',fp_out);
*/	
	}

	fputc('\n',fp_out);

	fclose(fp_out);
	fclose(fp_in);
}


void pure_composition()
{
	FILE *fp_in = fopen(infilename,"r");
	FILE *fp_out = fopen(outfilecomposition,"w");
	char ch;
	int i,len;

	int count[20];
	for(i=0; i<20; i++)
		count[i]=0;

	while((ch=fgetc(fp_in))!=EOF)
	{
		//the line containing > sign is rejected
		if(ch=='>')
		{
			while((ch=fgetc(fp_in))!='\n')
				;
			ch=fgetc(fp_in);
		}
	
		//initialize 20 amino acid counter for this line
	
		//ch = toupper(ch);
		//for the actual sequence, we compute the composition
		if(ch=='A')
			count[0]++;
		else if(ch=='C')
			count[1]++;
		else if(ch=='D')
			count[2]++;
		else if(ch=='E')
			count[3]++;
		else if(ch=='F')
			count[4]++;
		else if(ch=='G')
			count[5]++;
		else if(ch=='H')
			count[6]++;
		else if(ch=='I')
			count[7]++;
		else if(ch=='K')
			count[8]++;
		else if(ch=='L')
			count[9]++;
		else if(ch=='M')
			count[10]++;
		else if(ch=='N')
			count[11]++;
		else if(ch=='P')
			count[12]++;
		else if(ch=='Q')
			count[13]++;
		else if(ch=='R')
			count[14]++;
		else if(ch=='S')
			count[15]++;
		else if(ch=='T')
			count[16]++;
		else if(ch=='V')
			count[17]++;
		else if(ch=='W')
			count[18]++;
		else if(ch=='Y')
			count[19]++;

		if(ch=='\n')
		{
			//compute composition and write to file
			len=0;
			for(i=0; i<20; i++)
				len += count[i];
			for(i=0; i<20; i++)
				fprintf(fp_out,"%.1f\t",((float)(count[i]*100)/len));
			fprintf(fp_out,"\n");

			for(i=0; i<20; i++)
				count[i]=0;
		}
	}

	fclose(fp_out);
	fclose(fp_in);

}











